function [G,outMeanMom] = G_fun( Xi,Rho,Nu,M,dim )

% Produce moments from model at (Xi,Rho,Nu)
sim_moms            = produce_moments(Xi,Rho,Nu);

% Compute function G
G                   = M - mean(sim_moms.M,2);

if nargin > 4
    G = G(dim);
end

if nargout > 1
    outMeanMom = mean(sim_moms.M,2);
end
    
end

